library(psych)
library(ggpubr)
library(tidyverse)
library(lmerTest)

Introduction

Regions of Interest

King 2019 Atlas

We will use regions 1 and 2 as control regions since there are related to motor processing. For ROIs we will use regions 5, 6, 7, 8, 9, and 10 since they all overall with crus I and II which have been the main regions implicated in social cognition. There are also some studies which point to Lobule VI and Crus II’s involvement in reward processing.

Import Data

participants <- read.table(file='../participants_good.tsv', sep = '\t', header = TRUE)
participants

Import first level ROI data

all_sub_roi_con <- read.csv('../univariate_roi/all_sub_roi_contrasts.csv')
all_sub_roi_con$subject_id <- all_sub_roi_con$subj

# Use "adult as group label
all_sub_roi_con$group[all_sub_roi_con$group == 'control'] <- 'college'

all_sub_roi_con

Define hypothesis ROIs

hyp_rois <- c('TPJ_b-social','ATL_b-social','AMG_r-social','PCC_b-social','dmPFC_b-social',
              'vmPFC_b-reward',
              'region01','region02','region05','region06','region07',
              'region08','region09','region10', 
              'striatum_ventral','striatum_dorsala','striatum_dorsalp')

hyp_rois_names <- c('TPJ','ATL','AMG (r)','PCC','dmPFC',
                    'vmPFC',
                    'CB Region 1','CB Region 2','CB Region 5',
                    'CB Region 6','CB Region 7','CB Region 8','CB Region 9','CB Region 10',
                    'Ventral Striatum','Anterior \n Dorsal Striatum','Posterior \n Dorsal Striatum')

# Turn ROI column into a factor for sorting during visualization
all_sub_roi_con$roi <- factor(all_sub_roi_con$roi, levels=hyp_rois)

Positive Wins vs Loses

Young Adults

# Filter data for group, contrast, and ROIs
adult_data <- all_sub_roi_con[all_sub_roi_con$group %in% 'college', ]
adult_data_pos <- adult_data[adult_data$contrast %in% 'positive_winVlos', ]
adult_data_pos <- adult_data_pos[adult_data_pos$roi %in% hyp_rois, ]

# Filter for task
adult_data_pos_mdoors <- adult_data_pos[adult_data_pos$task %in% 'mdoors', ]
adult_data_pos_social <- adult_data_pos[adult_data_pos$task %in% 'social', ]

** Function: Conduct One Sample T-Tests for ROIs**

# One sample t-test
mass_t_tests <- function(df) {
  # Conduct one-sample t-tests to examine if ROI contrast is greater than 0
  df.t_tests = df %>%
               group_by(roi) %>%
               summarise(Statistic = t.test(contrast_mean, mu = 0)$statistic,
                         P = t.test(contrast_mean, mu = 0)$p.value,
                         Sig = ifelse(P < 0.05, "*", ""),
                         df = t.test(contrast_mean, mu = 0)$parameter,
                         Estimate = t.test(contrast_mean, mu = 0)$estimate,
                         CI_l = t.test(contrast_mean, mu = 0)$conf.int[1],
                         CI_u = t.test(contrast_mean, mu = 0)$conf.int[2],
                         MaxWidth = max(contrast_mean))
  
  # Conduct multiple comparisons correction
  stats <- p.adjust(df.t_tests$P,method='fdr')
  P_adjust <- stats
  Sig_adjust = ifelse(P_adjust < 0.005, "***", ifelse(P_adjust < 0.01, "**", 
                                               ifelse(P_adjust < 0.05, "*", "")))
  df.t_tests$P_adjust <- P_adjust
  df.t_tests$Sig_adjust <- Sig_adjust
  
  return(df.t_tests)
}


# Two sample t-test
mass_twos_t_tests <- function(df) {
  # Conduct one-sample t-tests to examine if ROI contrast is greater than 0
  df.t_tests <- data.frame(roi = character(),
                           Statistic = numeric(),
                           P = numeric(),
                           Sig = character(),
                           df = numeric(),
                           Estimate = numeric(),
                           CI_l = numeric(),
                           CI_u = numeric(),
                           MaxWidth = numeric(),
                           P_adjust = numeric(),
                           Sig_adjust = character())


  for (n in 1:length(hyp_rois)) {
    temp_df <- df[df$roi == hyp_rois[n],]
    group1 <- temp_df[temp_df$group == 'kid',]$contrast_mean
    group2 <- temp_df[temp_df$group == 'college',]$contrast_mean
    temp_ttest <- t.test(group1, group2)
    
    df.t_tests[n,]$roi <- hyp_rois[n]
    df.t_tests[n,]$Statistic <- temp_ttest$statistic
    df.t_tests[n,]$P <- temp_ttest$p.value
    df.t_tests[n,]$Sig <- ifelse(df.t_tests[n,]$P < 0.05, "*", "")
    df.t_tests[n,]$df <- temp_ttest$parameter
    df.t_tests[n,]$Estimate <- abs(temp_ttest$estimate[1] - temp_ttest$estimate[2])
    df.t_tests[n,]$CI_l <- temp_ttest$conf.int[1]
    df.t_tests[n,]$CI_u <- temp_ttest$conf.int[2]
    df.t_tests[n,]$MaxWidth <- max(temp_df$contrast_mean)
  }
    
  # Conduct multiple comparisons correction
  stats <- p.adjust(df.t_tests$P, method='fdr')
  P_adjust <- stats
  Sig_adjust = ifelse(P_adjust < 0.005, "***", ifelse(P_adjust < 0.01, "**", 
                                               ifelse(P_adjust < 0.05, "*", "")))
  df.t_tests$P_adjust <- P_adjust
  df.t_tests$Sig_adjust <- Sig_adjust
  
  return(df.t_tests)
}

Monetary Task

# Conduct one-sample t-tests to examine if ROI contrast is greater than 0
adult_data_pos_mdoors.t_tests <- mass_t_tests(adult_data_pos_mdoors)

# Plot
p <- ggbarplot(adult_data_pos_mdoors, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se',
               position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = adult_data_pos_mdoors.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5)) + 
  xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names)

ggsave('../univariate_roi/adult_pos_winVlos_ttest_bar_mdoors.png')

Social Task

# Conduct one-sample t-tests to examine if ROI contrast is greater than 0
adult_data_pos_social.t_tests <- mass_t_tests(adult_data_pos_social)


# Plot
p <- ggbarplot(adult_data_pos_social, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se',
               position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.35), size = 6,
              data = adult_data_pos_social.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5)) + 
  xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names)

ggsave('../univariate_roi/adult_pos_winVlos_ttest_bar_social.png')

Adolescents

# Filter data
adole_data <- all_sub_roi_con[all_sub_roi_con$group %in% 'kid', ]
adole_data_pos <- adole_data[adole_data$contrast %in% 'positive_winVlos', ]
adole_data_pos <- adole_data_pos[adole_data_pos$roi %in% hyp_rois, ]

adole_data_pos_mdoors <- adole_data_pos[adole_data_pos$task %in% 'mdoors', ]
adole_data_pos_social <- adole_data_pos[adole_data_pos$task %in% 'social', ]

Monetary Task

# One sample t-tests
adole_data_pos_mdoors.t_tests <- mass_t_tests(adole_data_pos_mdoors)


# Plot
p <- ggbarplot(adole_data_pos_mdoors, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se',
               position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = adole_data_pos_mdoors.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5)) + 
  xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names)

ggsave('../univariate_roi/adole_pos_winVlos_ttest_bar_mdoors.png')

Social Task

# One sample t-tests
adole_data_pos_social.t_tests <- mass_t_tests(adole_data_pos_social)


# Plot
p <- ggbarplot(adole_data_pos_social, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se',
               position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = adole_data_pos_social.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5)) + 
  xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names)

ggsave('../univariate_roi/adole_pos_winVlos_ttest_bar_social.png')

Group Differences

# Filter data
allg_data_pos <- all_sub_roi_con[all_sub_roi_con$contrast %in% 'positive_winVlos', ]
allg_data_pos <- allg_data_pos[allg_data_pos$roi %in% hyp_rois, ]

Monetary Task

allg_data_pos_mdoors <- allg_data_pos[allg_data_pos$task %in% 'mdoors', ]
# Conduct ANOVA
model_allg_pos_mdoors <- aov(contrast_mean ~ group*roi,
                           data = allg_data_pos_mdoors)
summary(model_allg_pos_mdoors)
##               Df Sum Sq Mean Sq F value Pr(>F)
## group          1   0.01 0.00743   0.097  0.755
## roi           16   1.56 0.09771   1.280  0.202
## group:roi     16   1.16 0.07246   0.949  0.512
## Residuals   1003  76.59 0.07636
s1 <- adult_data_pos_mdoors.t_tests$Sig_adjust
s2 <- adole_data_pos_mdoors.t_tests$Sig_adjust

for (n in 1:length(s1)) {
  if (s1[n] == "") {
    s1[n] = "#FFFFFF"
  } else {
    s1[n] = "#0072B2"
  }
  if (s2[n] == "") {
    s2[n] = "#FFFFFF"
  } else {
    s2[n] = "#D55E00"
  }
}

test <- c(rbind(s1, s2))
bar_fill_one_sample <- function(df, g1.t_tests, g2.t_tests) {
  df$plot_color <- ''
  for (i_roi in unique(df$roi)) {
    temp_sig <- g1.t_tests[g1.t_tests$roi == i_roi,'Sig_adjust']
    if (temp_sig == "") {
      temp_color <- 'adole_not'
    } else {
      temp_color <- 'adole_sig'
    }
    df[df$group == 'kid' & df$roi == i_roi,]$plot_color <- temp_color
    
    temp_sig <- g2.t_tests[g2.t_tests$roi == i_roi,'Sig_adjust']
    if (temp_sig == "") {
      temp_color <- 'adult_not'
    } else {
      temp_color <- 'adult_sig'
    }
    df[df$group == 'college' & df$roi == i_roi,]$plot_color <- temp_color
  }
  
  return(df)
} 
allg_data_pos_mdoors$plot_color <- ''
for (i_roi in unique(allg_data_pos_mdoors$roi)) {
  temp_sig <- adole_data_pos_mdoors.t_tests[adole_data_pos_mdoors.t_tests$roi == i_roi,'Sig_adjust']
  if (temp_sig == "") {
    temp_color <- 'adole_not'
  } else {
    temp_color <- 'adole_sig'
  }
  allg_data_pos_mdoors[allg_data_pos_mdoors$group == 'kid' & allg_data_pos_mdoors$roi == i_roi,]$plot_color <- temp_color
  
  temp_sig <- adult_data_pos_mdoors.t_tests[adult_data_pos_mdoors.t_tests$roi == i_roi,'Sig_adjust']
  if (temp_sig == "") {
    temp_color <- 'adult_not'
  } else {
    temp_color <- 'adult_sig'
  }
  allg_data_pos_mdoors[allg_data_pos_mdoors$group == 'college' & allg_data_pos_mdoors$roi == i_roi,]$plot_color <- temp_color
}
allg_data_pos_mdoors.t_tests <- mass_twos_t_tests(allg_data_pos_mdoors)

fill.colors <- c(adole_not = "#FFFFFF", adole_sig = "#7fb8d8", adult_not ="#FFFFFF", adult_sig = "#eaae7f")
outline.colors <- c("#D55E00", "#0072B2")

# Plot
p <- ggbarplot(allg_data_pos_mdoors, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se', add.params = list(group = "group"),
               group = 'group', color='group', position = position_dodge(0.8), 
               palette = outline.colors, fill = 'plot_color')
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = allg_data_pos_mdoors.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5), legend.position='right') + 
  xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names) + 
  scale_fill_manual(values=fill.colors) + 
  guides(fill = FALSE) + 
  scale_colour_manual(name = 'Group', values = outline.colors, 
                      labels=c("Young Adult", "Adolescent")) + 
  guides(color = guide_legend(override.aes = list(fill = outline.colors)))

ggsave('../univariate_roi/allp_pos_winVlos_anova_bar_mdoors.png')

Social Task

# Filter data
allg_data_pos_social <- allg_data_pos[allg_data_pos$task %in% 'social', ]
# Conduct ANOVA
model_allg_pos_social <- aov(contrast_mean ~ group*roi,
                           data = allg_data_pos_social)
summary(model_allg_pos_social)
##               Df Sum Sq Mean Sq F value  Pr(>F)   
## group          1   0.50  0.5035   7.921 0.00498 **
## roi           16   1.89  0.1183   1.862 0.02045 * 
## group:roi     16   0.51  0.0320   0.503 0.94649   
## Residuals   1003  63.76  0.0636                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
allg_data_pos_social.t_tests <- mass_twos_t_tests(allg_data_pos_social)

allg_data_pos_social <- bar_fill_one_sample(allg_data_pos_social, 
                                            adole_data_pos_social.t_tests,
                                            adult_data_pos_social.t_tests)

fill.colors <- c(adole_not = "#FFFFFF", adole_sig = "#7fb8d8", adult_not ="#FFFFFF", adult_sig = "#eaae7f")
outline.colors <- c("#D55E00", "#0072B2")

# Plot
p <- ggbarplot(allg_data_pos_social, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se', add.params = list(group = "group"),
               group = 'group', color='group', position = position_dodge(0.8), 
               palette = outline.colors, fill = 'plot_color')
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = allg_data_pos_social.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5), legend.position='right') + 
  xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names) + 
  scale_fill_manual(values=fill.colors) + 
  guides(fill = FALSE) + 
  scale_colour_manual(name = 'Group', values = outline.colors, 
                      labels=c("Young Adult", "Adolescent")) + 
  guides(color = guide_legend(override.aes = list(fill = outline.colors)))

ggsave('../univariate_roi/allp_pos_winVlos_anova_bar_social.png')

All Wins vs Loses

Young Adults

# Filter data
adult_data <- all_sub_roi_con[all_sub_roi_con$group %in% 'college', ]
adult_data_all <- adult_data[adult_data$contrast %in% 'all_winVlos', ]
adult_data_all <- adult_data_all[adult_data_all$roi %in% hyp_rois, ]

adult_data_all_mdoors <- adult_data_all[adult_data_all$task %in% 'mdoors', ]
adult_data_all_social <- adult_data_all[adult_data_all$task %in% 'social', ]

Monetary Task

# Conduct t-tests
adult_data_all_mdoors.t_tests <- mass_t_tests(adult_data_all_mdoors)


# Plot
p <- ggbarplot(adult_data_all_mdoors, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se',
               position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = adult_data_all_mdoors.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5)) + 
  xlab("Regions of Interests") + ylab("Mean Response (all wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names)

ggsave('../univariate_roi/adult_all_winVlos_ttest_bar_mdoors.png')

Social Task

# T-test
adult_data_all_social.t_tests <- mass_t_tests(adult_data_all_social)


# Plot
p <- ggbarplot(adult_data_all_social, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se',
               position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = adult_data_all_social.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5)) + 
  xlab("Regions of Interests") + ylab("Mean Response (all wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names)

ggsave('../univariate_roi/adult_all_winVlos_ttest_bar_social.png')

Adolescents

# Filter data
adole_data <- all_sub_roi_con[all_sub_roi_con$group %in% 'kid', ]
adole_data_all <- adole_data[adole_data$contrast %in% 'all_winVlos', ]
adole_data_all <- adole_data_all[adole_data_all$roi %in% hyp_rois, ]

adole_data_all_mdoors <- adole_data_all[adole_data_all$task %in% 'mdoors', ]
adole_data_all_social <- adole_data_all[adole_data_all$task %in% 'social', ]

Monetary Task

# T-test
adole_data_all_mdoors.t_tests <- mass_t_tests(adole_data_all_mdoors)


# Plot
p <- ggbarplot(adole_data_all_mdoors, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se',
               position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = adole_data_all_mdoors.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5)) + 
  xlab("Regions of Interests") + ylab("Mean Response (all wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names)

ggsave('../univariate_roi/adole_all_winVlos_ttest_bar_mdoors.png')

Social Task

# T-test
adole_data_all_social.t_tests <- mass_t_tests(adole_data_all_social)


# Plot
p <- ggbarplot(adole_data_all_social, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se',
               position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = adole_data_all_social.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5)) + 
  xlab("Regions of Interests") + ylab("Mean Response (all wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names)

ggsave('../univariate_roi/adole_all_winVlos_ttest_bar_social.png')

Group Differences

# Filter data
allg_data_all <- all_sub_roi_con[all_sub_roi_con$contrast %in% 'all_winVlos', ]
allg_data_all <- allg_data_all[allg_data_all$roi %in% hyp_rois, ]

Monetary Task

allg_data_all_mdoors <- allg_data_all[allg_data_all$task %in% 'mdoors', ]
# ANOVA
model_allg_all_mdoors <- aov(contrast_mean ~ group*roi,
                           data = allg_data_all_mdoors)
summary(model_allg_all_mdoors)
##               Df Sum Sq Mean Sq F value Pr(>F)
## group          1   0.07 0.06990   0.815  0.367
## roi           16   0.85 0.05288   0.616  0.873
## group:roi     16   1.28 0.07986   0.931  0.533
## Residuals   1003  86.05 0.08580
allg_data_all_mdoors.t_tests <- mass_twos_t_tests(allg_data_all_mdoors)

allg_data_all_mdoors <- bar_fill_one_sample(allg_data_all_mdoors, 
                                            adole_data_all_mdoors.t_tests,
                                            adult_data_all_mdoors.t_tests)

fill.colors <- c(adole_not = "#FFFFFF", adole_sig = "#7fb8d8", adult_not ="#FFFFFF", adult_sig = "#eaae7f")
outline.colors <- c("#D55E00", "#0072B2")

# Plot
p <- ggbarplot(allg_data_all_mdoors, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se', add.params = list(group = "group"),
               group = 'group', color='group', position = position_dodge(0.8), 
               palette = outline.colors, fill = 'plot_color')
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = allg_data_all_mdoors.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5), legend.position='right') + 
  xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names) + 
  scale_fill_manual(values=fill.colors) + 
  guides(fill = FALSE) + 
  scale_colour_manual(name = 'Group', values = outline.colors, 
                      labels=c("Young Adult", "Adolescent")) + 
  guides(color = guide_legend(override.aes = list(fill = outline.colors)))

ggsave('../univariate_roi/allp_all_winVlos_anova_bar_mdoors.png')

Social Task

# Filter data
allg_data_all_social <- allg_data_all[allg_data_all$task %in% 'social', ]
# ANOVA
model_allg_all_social <- aov(contrast_mean ~ group*roi,
                           data = allg_data_all_social)
summary(model_allg_all_social)
##               Df Sum Sq Mean Sq F value   Pr(>F)    
## group          1   0.90  0.9011  14.057 0.000188 ***
## roi           16   1.76  0.1101   1.718 0.038252 *  
## group:roi     16   0.17  0.0109   0.170 0.999908    
## Residuals   1003  64.29  0.0641                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
allg_data_all_social.t_tests <- mass_twos_t_tests(allg_data_all_social)

allg_data_all_social <- bar_fill_one_sample(allg_data_all_social, 
                                            adole_data_all_social.t_tests,
                                            adult_data_all_social.t_tests)

fill.colors <- c(adole_not = "#FFFFFF", adole_sig = "#7fb8d8", adult_not ="#FFFFFF", adult_sig = "#eaae7f")
outline.colors <- c("#D55E00", "#0072B2")

# Plot
p <- ggbarplot(allg_data_all_social, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se', add.params = list(group = "group"),
               group = 'group', color='group', position = position_dodge(0.8), 
               palette = outline.colors, fill = 'plot_color')
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = allg_data_all_social.t_tests) + 
  theme(axis.text.x = element_text(angle =90, vjust = 0.5), legend.position='right') + 
  xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") + 
  scale_x_discrete(labels=hyp_rois_names) + 
  scale_fill_manual(values=fill.colors) + 
  guides(fill = FALSE) + 
  scale_colour_manual(name = 'Group', values = outline.colors, 
                      labels=c("Young Adult", "Adolescent")) + 
  guides(color = guide_legend(override.aes = list(fill = outline.colors)))

ggsave('../univariate_roi/allp_all_winVlos_anova_bar_social.png')

Combined Plots

# Plot
p1 <- ggbarplot(allg_data_pos_mdoors, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se', add.params = list(group = "group"),
               group = 'group', color='group', position = position_dodge(0.8), 
               palette = outline.colors, fill = 'plot_color') + 
  geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = allg_data_pos_mdoors.t_tests) + 
  #theme(axis.text.x = element_text(angle =90, vjust = 0.5), legend.position='right') + 
  theme(axis.text.x=element_blank(), axis.title.x = element_blank(),
        plot.title = element_text(face='bold')) + 
  #xlab("\n") + 
  ylab("Mean Response") + 
  ggtitle('A. Positive Wins > Positive Losses Contrast') + 
  #scale_x_discrete(labels=hyp_rois_names) + 
  scale_fill_manual(values=fill.colors) + 
  guides(fill = FALSE) + 
  scale_colour_manual(name = 'Group', values = outline.colors, 
                      labels=c("Young Adult", "Adolescent")) + 
  guides(color = guide_legend(override.aes = list(fill = outline.colors)))


p2 <- ggbarplot(allg_data_all_mdoors, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se', add.params = list(group = "group"),
               group = 'group', color='group', position = position_dodge(0.8), 
               palette = outline.colors, fill = 'plot_color') +
  geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = allg_data_all_mdoors.t_tests) + 
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5), legend.position='right',
        plot.title = element_text(face='bold')) + 
  xlab("Regions of Interests") + ylab("Mean Response") + 
  ggtitle('B. All Wins > All Losses Contrast') + 
  scale_x_discrete(labels=hyp_rois_names) + 
  scale_fill_manual(values=fill.colors) + 
  guides(fill = FALSE) + 
  scale_colour_manual(name = 'Group', values = outline.colors, 
                      labels=c("Young Adult", "Adolescent")) + 
  guides(color = guide_legend(override.aes = list(fill = outline.colors)))


multi_plot <- ggarrange(p1, p2, nrow = 2, common.legend = TRUE, legend = 'right',
          heights = c(5,6.75), align = 'v') 

annotate_figure(multi_plot, top = text_grob("ROI Responses for Monetary Task\n", 
               face = "bold", size = 18))

ggsave('../univariate_roi/allp_anova_bar_mdoors.png')
# Plot
p1 <- ggbarplot(allg_data_pos_social, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se', add.params = list(group = "group"),
               group = 'group', color='group', position = position_dodge(0.8), 
               palette = outline.colors, fill = 'plot_color') + 
  geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = allg_data_pos_social.t_tests) + 
  #theme(axis.text.x = element_text(angle =90, vjust = 0.5), legend.position='right') + 
  theme(axis.text.x=element_blank(), axis.title.x = element_blank(),
        plot.title = element_text(face='bold')) + 
  #xlab("\n") + 
  ylab("Mean Response") + 
  ggtitle('A. Positive Wins > All Losses Contrast') + 
  #scale_x_discrete(labels=hyp_rois_names) + 
  scale_fill_manual(values=fill.colors) + 
  guides(fill = FALSE) + 
  scale_colour_manual(name = 'Group', values = outline.colors, 
                      labels=c("Young Adult", "Adolescent")) + 
  guides(color = guide_legend(override.aes = list(fill = outline.colors)))


p2 <- ggbarplot(allg_data_all_social, x = 'roi', 
               y = 'contrast_mean', add = 'mean_se', add.params = list(group = "group"),
               group = 'group', color='group', position = position_dodge(0.8), 
               palette = outline.colors, fill = 'plot_color') +
  geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
              data = allg_data_all_social.t_tests) + 
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5), legend.position='right',
        plot.title = element_text(face='bold')) + 
  xlab("Regions of Interests") + ylab("Mean Response") + 
  ggtitle('B. All Wins > All Losses Contrast') + 
  scale_x_discrete(labels=hyp_rois_names) + 
  scale_fill_manual(values=fill.colors) + 
  guides(fill = FALSE) + 
  scale_colour_manual(name = 'Group', values = outline.colors, 
                      labels=c("Young Adult", "Adolescent")) + 
  guides(color = guide_legend(override.aes = list(fill = outline.colors)))


multi_plot <- ggarrange(p1, p2, nrow = 2, common.legend = TRUE, legend = 'right',
          heights = c(5,6.75), align = 'v') 

annotate_figure(multi_plot, top = text_grob("ROI Responses for Social Task\n", 
               face = "bold", size = 18))

ggsave('../univariate_roi/allp_anova_bar_social.png')

Stats Table

Create a stats table to summarize all univariate ROI analyses

Combine one-sample t-test data

adult_data_pos_mdoors.t_tests$group <- 'adult'
adult_data_pos_mdoors.t_tests$task <- 'mdoors'
adult_data_pos_mdoors.t_tests$contrast <- 'positive'

adult_data_pos_social.t_tests$group <- 'adult'
adult_data_pos_social.t_tests$task <- 'social'
adult_data_pos_social.t_tests$contrast <- 'positive'

adult_data_all_mdoors.t_tests$group <- 'adult'
adult_data_all_mdoors.t_tests$task <- 'mdoors'
adult_data_all_mdoors.t_tests$contrast <- 'all'

adult_data_all_social.t_tests$group <- 'adult'
adult_data_all_social.t_tests$task <- 'social'
adult_data_all_social.t_tests$contrast <- 'all'


adole_data_pos_mdoors.t_tests$group <- 'adolescent'
adole_data_pos_mdoors.t_tests$task <- 'mdoors'
adole_data_pos_mdoors.t_tests$contrast <- 'positive'

adole_data_pos_social.t_tests$group <- 'adolescent'
adole_data_pos_social.t_tests$task <- 'social'
adole_data_pos_social.t_tests$contrast <- 'positive'

adole_data_all_mdoors.t_tests$group <- 'adolescent'
adole_data_all_mdoors.t_tests$task <- 'mdoors'
adole_data_all_mdoors.t_tests$contrast <- 'all'

adole_data_all_social.t_tests$group <- 'adolescent'
adole_data_all_social.t_tests$task <- 'social'
adole_data_all_social.t_tests$contrast <- 'all'


combined_ttests <- rbind(#adult_data_pos_mdoors.t_tests,
                         #adult_data_pos_social.t_tests,
                         adult_data_all_mdoors.t_tests,
                         adult_data_all_social.t_tests,
                         #adole_data_pos_mdoors.t_tests,
                         #adole_data_pos_social.t_tests, 
                         adole_data_all_mdoors.t_tests,
                         adole_data_all_social.t_tests)

Fix up dataframe for pretty table

# Sort dataframe
one_sample_table <- combined_ttests[order(combined_ttests$group, 
                                          combined_ttests$task, 
                                          combined_ttests$roi), ]

# Rename ROI names
rep_str = c('TPJ_b-social'='TPJ', 'ATL_b-social'='ATL', 'AMG_r-social'='AMG',
            'PCC_b-social'='PCC', 'dmPFC_b-social'='dmPFC', 
            'vmPFC_b-reward'='vmPFC', 'region01'='CB Region 1',
            'region02'='CB Region 2', 'region05'='CB Region 5', 
            'region06'='CB Region 6', 'region07'='CB Region 7', 
            'region08'='CB Region 8', 'region09'='CB Region 9',
            'region10'='CB Region 10', 
            'striatum_dorsalp'='Posterior Dorsal Striatum', 
            'striatum_dorsala'='Posterior Anterior Striatum', 
            'striatum_ventral'='Ventral Striatum')
one_sample_table$roi <- str_replace_all(one_sample_table$roi, rep_str)

# Turn columns numeric
one_sample_table$Statistic <- as.numeric(one_sample_table$Statistic)
one_sample_table$Estimate <- as.numeric(one_sample_table$Estimate)
one_sample_table$CI_l <- as.numeric(one_sample_table$CI_l)
one_sample_table$CI_u <- as.numeric(one_sample_table$CI_u)
one_sample_table$P_adjust <- as.numeric(one_sample_table$P_adjust)

# Round numbers
one_sample_table$Statistic <- round(one_sample_table$Statistic,2)
one_sample_table$Estimate <- round(one_sample_table$Estimate,2)
one_sample_table$CI_l <- round(one_sample_table$CI_l,2)
one_sample_table$CI_u <- round(one_sample_table$CI_u,2)
one_sample_table$P_adjust <- round(one_sample_table$P_adjust,3)

# Combine confidence inverval ranges
one_sample_table$CI <- paste('[', as.character(one_sample_table$CI_l),
                              ', ', as.character(one_sample_table$CI_u), ']')

# Reorder columns
one_sample_table <- one_sample_table[, c('group', 'task', "roi", "Estimate", 
                                         "df", "Statistic",
                                           "P_adjust", "CI")]

colnames(one_sample_table) <- c('Group', 'Task', 'ROI','Mean','df','t-statistic','p','CI')
one_sample_table %>%
  kbl(caption = "RSA results for a full and partial similarity analysis") %>%
  kable_classic(full_width = F, html_font = "Times New Roman") %>%
  save_kable("mass_one_sample_ttests.html")


one_sample_table %>%
  kbl(caption = "RSA results for a full and partial similarity analysis") %>%
  kable_classic(full_width = F, html_font = "Times New Roman")
RSA results for a full and partial similarity analysis
Group Task ROI Mean df t-statistic p CI
adolescent mdoors TPJ 0.04 31 0.73 0.731 [ -0.07 , 0.15 ]
adolescent mdoors ATL 0.04 31 0.90 0.713 [ -0.05 , 0.12 ]
adolescent mdoors AMG 0.00 31 0.07 0.990 [ -0.1 , 0.1 ]
adolescent mdoors PCC 0.00 31 0.02 0.990 [ -0.19 , 0.2 ]
adolescent mdoors dmPFC 0.00 31 0.01 0.990 [ -0.12 , 0.12 ]
adolescent mdoors vmPFC 0.09 31 1.16 0.605 [ -0.07 , 0.26 ]
adolescent mdoors CB Region 1 0.01 31 0.23 0.990 [ -0.08 , 0.1 ]
adolescent mdoors CB Region 2 -0.01 31 -0.22 0.990 [ -0.11 , 0.09 ]
adolescent mdoors CB Region 5 0.09 31 1.64 0.605 [ -0.02 , 0.2 ]
adolescent mdoors CB Region 6 0.06 31 1.32 0.605 [ -0.03 , 0.15 ]
adolescent mdoors CB Region 7 0.07 31 1.20 0.605 [ -0.05 , 0.19 ]
adolescent mdoors CB Region 8 0.08 31 1.31 0.605 [ -0.05 , 0.21 ]
adolescent mdoors CB Region 9 0.09 31 1.54 0.605 [ -0.03 , 0.2 ]
adolescent mdoors CB Region 10 0.04 31 1.09 0.605 [ -0.04 , 0.12 ]
adolescent mdoors Ventral Striatum 0.09 31 1.77 0.605 [ -0.01 , 0.19 ]
adolescent mdoors Posterior Anterior Striatum 0.05 31 0.73 0.731 [ -0.08 , 0.17 ]
adolescent mdoors Posterior Dorsal Striatum 0.02 31 0.33 0.990 [ -0.09 , 0.13 ]
adolescent social TPJ 0.07 31 1.92 0.116 [ 0 , 0.14 ]
adolescent social ATL -0.04 31 -1.55 0.150 [ -0.1 , 0.01 ]
adolescent social AMG 0.01 31 0.25 0.808 [ -0.07 , 0.09 ]
adolescent social PCC 0.13 31 2.07 0.108 [ 0 , 0.26 ]
adolescent social dmPFC 0.09 31 1.75 0.138 [ -0.01 , 0.19 ]
adolescent social vmPFC 0.11 31 2.03 0.108 [ 0 , 0.22 ]
adolescent social CB Region 1 0.05 31 1.41 0.179 [ -0.02 , 0.12 ]
adolescent social CB Region 2 0.06 31 2.22 0.108 [ 0.01 , 0.12 ]
adolescent social CB Region 5 0.06 31 1.64 0.150 [ -0.02 , 0.14 ]
adolescent social CB Region 6 0.06 31 2.04 0.108 [ 0 , 0.13 ]
adolescent social CB Region 7 0.09 31 1.89 0.116 [ -0.01 , 0.18 ]
adolescent social CB Region 8 0.14 31 3.33 0.019 [ 0.06 , 0.23 ]
adolescent social CB Region 9 0.13 31 3.48 0.019 [ 0.05 , 0.2 ]
adolescent social CB Region 10 0.06 31 2.48 0.080 [ 0.01 , 0.12 ]
adolescent social Ventral Striatum 0.14 31 2.91 0.037 [ 0.04 , 0.25 ]
adolescent social Posterior Anterior Striatum 0.09 31 1.56 0.150 [ -0.03 , 0.22 ]
adolescent social Posterior Dorsal Striatum 0.07 31 1.54 0.150 [ -0.02 , 0.17 ]
adult mdoors TPJ 0.06 28 1.47 0.515 [ -0.02 , 0.14 ]
adult mdoors ATL 0.02 28 0.43 0.712 [ -0.06 , 0.09 ]
adult mdoors AMG 0.04 28 1.06 0.632 [ -0.04 , 0.12 ]
adult mdoors PCC 0.18 28 2.04 0.333 [ 0 , 0.36 ]
adult mdoors dmPFC 0.05 28 0.85 0.632 [ -0.07 , 0.17 ]
adult mdoors vmPFC 0.05 28 0.89 0.632 [ -0.07 , 0.17 ]
adult mdoors CB Region 1 -0.02 28 -0.49 0.712 [ -0.09 , 0.05 ]
adult mdoors CB Region 2 -0.02 28 -0.49 0.712 [ -0.1 , 0.06 ]
adult mdoors CB Region 5 -0.03 28 -0.77 0.633 [ -0.12 , 0.06 ]
adult mdoors CB Region 6 -0.04 28 -1.13 0.632 [ -0.11 , 0.03 ]
adult mdoors CB Region 7 -0.02 28 -0.46 0.712 [ -0.12 , 0.08 ]
adult mdoors CB Region 8 0.05 28 0.84 0.632 [ -0.07 , 0.17 ]
adult mdoors CB Region 9 -0.01 28 -0.15 0.881 [ -0.12 , 0.1 ]
adult mdoors CB Region 10 -0.03 28 -0.93 0.632 [ -0.1 , 0.04 ]
adult mdoors Ventral Striatum 0.08 28 2.06 0.333 [ 0 , 0.16 ]
adult mdoors Posterior Anterior Striatum 0.07 28 1.83 0.333 [ -0.01 , 0.15 ]
adult mdoors Posterior Dorsal Striatum 0.05 28 1.91 0.333 [ 0 , 0.1 ]
adult social TPJ 0.13 28 2.77 0.015 [ 0.03 , 0.23 ]
adult social ATL 0.10 28 2.85 0.014 [ 0.03 , 0.17 ]
adult social AMG 0.10 28 1.93 0.064 [ -0.01 , 0.2 ]
adult social PCC 0.19 28 2.00 0.058 [ 0 , 0.38 ]
adult social dmPFC 0.15 28 2.39 0.029 [ 0.02 , 0.28 ]
adult social vmPFC 0.15 28 2.15 0.046 [ 0.01 , 0.29 ]
adult social CB Region 1 0.09 28 2.67 0.016 [ 0.02 , 0.16 ]
adult social CB Region 2 0.08 28 2.71 0.016 [ 0.02 , 0.14 ]
adult social CB Region 5 0.13 28 3.21 0.008 [ 0.05 , 0.21 ]
adult social CB Region 6 0.12 28 3.57 0.004 [ 0.05 , 0.2 ]
adult social CB Region 7 0.14 28 3.11 0.009 [ 0.05 , 0.23 ]
adult social CB Region 8 0.21 28 4.61 0.001 [ 0.12 , 0.31 ]
adult social CB Region 9 0.19 28 4.33 0.001 [ 0.1 , 0.28 ]
adult social CB Region 10 0.10 28 2.84 0.014 [ 0.03 , 0.17 ]
adult social Ventral Striatum 0.19 28 4.94 0.001 [ 0.11 , 0.26 ]
adult social Posterior Anterior Striatum 0.16 28 3.58 0.004 [ 0.07 , 0.25 ]
adult social Posterior Dorsal Striatum 0.12 28 3.80 0.003 [ 0.06 , 0.19 ]

Supplementary Analyses

Group x Task x ROI Analysis

Positive Wins vs Losses

# Filter data
allg_data_pos <- all_sub_roi_con[all_sub_roi_con$contrast %in% 'positive_winVlos', ]
allg_data_pos <- allg_data_pos[allg_data_pos$roi %in% hyp_rois, ]
# ANOVA
model_allg_pos <- aov(contrast_mean ~ group*task*roi,
                           data = allg_data_pos)
summary(model_allg_pos)
##                  Df Sum Sq Mean Sq F value   Pr(>F)    
## group             1   0.32  0.3166   4.526 0.033509 *  
## task              1   0.65  0.6497   9.286 0.002339 ** 
## roi              16   2.94  0.1837   2.626 0.000433 ***
## group:task        1   0.19  0.1943   2.777 0.095778 .  
## group:roi        16   0.76  0.0476   0.681 0.815591    
## task:roi         16   0.52  0.0323   0.462 0.964770    
## group:task:roi   16   0.91  0.0568   0.812 0.672752    
## Residuals      2006 140.35  0.0700                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Plot
p <- ggplot(data = allg_data_pos, aes(x=roi, y=contrast_mean, fill=task, alpha=group))
p + geom_bar(stat = 'summary', fun = 'mean', position=position_dodge()) +
  scale_alpha_manual(values = c(0.3, 1)) + 
  theme_classic() + 
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5))

ggsave('../univariate_roi/allp_pos_winVlos_anova_bar.png')

All Wins vs Losses

# Filter
allg_data_all <- all_sub_roi_con[all_sub_roi_con$contrast %in% 'all_winVlos', ]
allg_data_all <- allg_data_all[allg_data_all$roi %in% hyp_rois, ]
# ANOVA
model_allg_all <- aov(contrast_mean ~ group*task*roi,
                           data = allg_data_all)
summary(model_allg_all)
##                  Df Sum Sq Mean Sq F value   Pr(>F)    
## group             1   0.23  0.2345   3.129  0.07706 .  
## task              1   2.57  2.5675  34.257 5.63e-09 ***
## roi              16   2.23  0.1393   1.858  0.02010 *  
## group:task        1   0.74  0.7364   9.826  0.00175 ** 
## group:roi        16   0.76  0.0473   0.631  0.86084    
## task:roi         16   0.38  0.0237   0.316  0.99538    
## group:task:roi   16   0.70  0.0435   0.580  0.90098    
## Residuals      2006 150.35  0.0749                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
p <- ggplot(data = allg_data_all, aes(x=roi, y=contrast_mean, fill=task, alpha=group))
p + geom_bar(stat = 'summary', fun = 'mean', position=position_dodge()) +
  scale_alpha_manual(values = c(0.3, 1)) + 
  theme_classic() + 
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5))

ggsave('../univariate_roi/allp_all_winVlos_anova_bar.png')